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Kecurslvo  aleorithns  for  the  solution  of  lin- 
Wfc least-squares  estimation  problems  liavc  boon 
basoil  mainly  on  .sta tc-spaco  nodols.  It  lias  boon 
knov.n,  hov.ovor,  that  such  aleorithns  exist  for 
stationary  tiino-scr.Los , using  input-output  de- 
scriptions (c.g.,  covariance  matrices).  Wo  in- 
troduco  a way  of  classifying  stochastic  processes 
la  terms  of  their  "distance"  from  stationarity 
that  loads  to  a derivation  of  an  efficient  Lcvln- 
son-typo  algorithm  for  arbitrary  (nonstationary) 
processes.  By  adding  structure  to  the  covariance 
ii:atrix,  these  general  results  specialize  to  state- 
space  typo  estimation  algorithms.  I.a  particular, 
the  Chandraselduir  equations  arc  shoivn  to  be  the 
natural  descendants  of  the  Ixivinson  algorithm. 


1.  Introduction 

The  problem  of  linear  least  squares  estima- 
tion has  been  studied  extensively  and  various 
ce'hotls  of  solution  have  been  developed.  Those 
may  bo  classified  into  estimation  algorithms  de- 
rived from  input-output  data  or  from  other  "ex- 
ternal" syatem  descriptions  and  algorithms  derived 
from  state-space  or  "internal"  models.  In  the  last 
decade  the  field  of  linear  least-squares  estima- 
tion has  been  dominated  by  state-estimation,  in 
particular  by  the  recursive  Kalman-Bucy  filter 
algorithm  and  its  various  versions,  v.hich  rely 
heavily  on  the  availability  of  state-space  models. 

In  many  applications,  however,  a state-space 
model  is  not  readily  available,  and  it  would  be 
preferable  to  have  algorithms  that  use  dii’cctly 
the  covariance  information  of  the  observed  pro- 
cess. The  solution  of  the  estimation  problem  is 
closely  related  to  the  pi-oblcm  of  inverting  the 
covariance  matrix.  Therefore,  tho  computational 
efficiency  of  escimatlon  algoritlims  i.:  strongly 
dependent  on  the  amount  of  computation  required 
for  inverting  an  appropriate  matrix.  For  illu- 
"^tralion  we  .sluill  mention  tho  i.mportant  cxanplo 
f”of  a stationary  process  and  Its  'rooplltz-t>pc 
^ “t.ovar  isTnee  matrix.  If  has  hcen  shown  [1-3]  that 
making  use  of  tlic  S[K'cial  structure  of  .au  N XN 
l-5roepl  1 1-’.  matrix,  It  can  be  Invcrtid  in  0(N~)  op- 
cr.itions  (null  I I’l  iuat  foils  and  additions),  compared 
lUo  0(.N'*)  oper.il  ions  required  In  ipniocal  for  llie 
— 4nvc'r:  if>ii  of  an  arbitrary  imilrlx.  It  lias  also 


been  long  known  In  certain  fields  (c.g.,  geophys- 
ical data  pi’ocessing  [1]  and  speech  compression 
studies  [5])  that  recursive  solutions  can  bo  ob- 
tained for  tho  prediction  of  stationary  processes. 
In  particular,  tlio  so-called  Ix-vinson  algorithm 
computes  tho  optimal-predictor  In  0(N“)  operations. 

It  has  seemed  In  the  past  tliat  tho  prediction 
of  r.onstat ionary  processes  would  require  0(N'3)  op- 
ci-ations  unless  wo  can  impose  a state-space  struc- 
ture on  tho  signal .and  noiso  processes.  However, 

It  is  not  unreasonable  to  e.xpect  that  between  tho 
highly  structured  Tocplltz  matrix  (or  stationary 
process)  and  a completely  arbitrary  covai-ianco  ma- 
trt.x,  th.ore  sliould  exist  matrices  (or  pcocesscs) 
witli  varying  degrees  of  .structure  and  tliat  tiiis 
structure  could  bo  somehow  utilized  In  reducing 
tho  amount  of  computation  involved  in  the  estima- 
tion problem.  That  this  Is  indeed  possible  has 
been  first  shown  in  [6]  by  introducing  the  concept 
of  "shift  (low)  rank,"  see  also  [7-9],  and  subse- 
quently in  [10,11].  These  results  wore  motivated 
for  discrete-time  problems  by  the  work  of  Levinson 
[14]  and  Golub  [IS],  and  for  continuous-time  by 
tho  Cliar.drascldiar-typo  equations  and  their  further 
developments  in  [12,13], 

In  this  papor  we  shall  introduce  Q,  an  index 
of  "distance  from  stationarity"  of  an  arbitrary 
nonstationary  process.  We  shall  show  how  recur- 
sive solutions  requiring  of  order  ON“  operations 
can  bo  obtained  for  such  processes  with  or  without 
assuming  a state-space  structure.  In  tiic  station- 
ary case  our  solution  reduces  to  known  algorith.rs 
given  in  [14,15]. 

Finally,  we  shall  show  tliat,  if  the  not  nec- 
essarily stationary  processes  arc  known  to  cone 
from  stato-space  m.odcls,  tiion  this  additional 
structural  information  can  be  used  to  reduce  our 
general  solution  algorithm  to  tlio  previously  knowr. 
Cii.indrasokfi.'ir-typo  oquation.s.  Ttiis  means  tli.'.t  we 
liavo  been  able  to  properly  imbed  the  state-space 
assumption  into  a general  input-output  framework. 


2.  A General  Linear  F,;.tir.ation  Problem 

V,e  shall  con.ilder  the  problem  of  cstiiratlng 
stocliartlc  process  .\(0  from  observation:  of 

related  pri>ee.'..s  >(•).  1/>L  >(■)  ( p-d  ini-n.;  lor..''.  I ) 

and  .x(-)  (n-sl  It  lens  loii.i  1 ) have  covarl.ince  r.itii- 
ces 
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t 

or  In  matrix  form 


0 < 1,J  < N ; 
0 < l.J  < N ; 


The  beski^near  least  squares  estimate  of  x^j 
given  p^t  - 0 < 1 < N - 1 ) has  the  form* 


N-1 


x(N  N - 1) 


1=0 


h (N,l)  y, 
xy  1 


The  optimal  one-step  ahead  predictor  hjj„(N,')  can 
be  determined  by  using  the  well  known  orthogonal- 
ity condition  on  the  prediction  error 


Xj^,  - x(n]n  - 1)  j.  y^^  ; 


which  means 


0 < k < N - 1 


N-1 


0 = e[(x^,-?(nIn-d)  *’i.v 


or  In  matrix  form 


- Ko 

where  h”  A [h  (N,0)  , . . . ,hjjy(N,N  - 1)  ] , an  n XNp 


matrix. 

Note  that  by  setting  x = y,  we  get  an  equa- 
tion defining  the  predictor  h(*,’)  of  the  ob- 
served process  itself,  i.e., 


therefore 


?(n|x  - 1)  = ^ h(N,i)  y . 
1=0 

“ [^N.O'  '"n.N-i]  ■ 

This  can  bo  rewritten  as 

= [0 0,E]  , 


(2) 


(3) 


where  E Is  determined  by  the  left-hand  side  of 

this  cq\iatlon.  For  estimating  x(0  at  a time 

Instant  within  the  observation  Interval  (0,N),  wo 

have  to  find  the  optlnml  filter  (smoother) 

II  ( • , • :N) , where 
xy 

N, 

x(k|x)  = ^ H (k,l;N)  y . 

1=0 

Using  the  orthogonality  condition  again,  nowon  the 
form 


X|^  - x(kjN)  j.  y^  , 


0 < 1 < N 


we  got 
N 


y II  (k,l  ;N)  r.  ^ , 0 < kf  < V 

1,0  k,o  ' - - 


I.^*  R^  = R^ 
xy  xy 


<4) 


H**  A [ii  (i,J;N)1  , 
xy  = L xy  ’ J 


0 < l,j  < N . 


Tlie  last  two  matrix  equations  Illustrate  the  fact 
that  solving  the  estimation  problem  Is  closely  re- 
lated to  the  problem  of  Inverting  the  covariance 
matrix,  since  for  both  solutions  of  the  predic- 
tion- and  the  smoothing-problem  we  got 

- [“"f ' fc "> 

= r"  [r"!"'  . 
xy  xy  L J 


(6) 


3.  Structural  Assumptions 

To  aid  In  specifying  the  structural  assump- 
tions, we  shall  introduce  some  notation  for  the 
"shifted -difference"  operators  6[*]  and  J[*l 
that  play  a central  role  in  this  paper. 


Let, 


^ K.j] ' 


0 < i,J  < N , 


where  are  pXp  matrices  (i.e.,  S is  a 

block  matrix).  Then  define 


B[S]  0 


J[S) 


®1,1  •••  ®1,N 

®0,0  ®0,N-1 

• 

_ 

• 

• • 

• • 

_®N,1  ®N,N_ 

_°N-1,0  ®N-1,N-1_ 

®0,0  •••  ®0,N 

“ 01 

• 

I ®o,o  “o. 

* 1 

. 

* • 

®N,0  •••  ®N,N 

0i 

1 ®N-1,0  •••  ®N-1 

N,0 


_ 2>!5 

B[S] 


Wo  can  now  define  the  (block)  displacement  rank  Q 
of  matrix  S ns 


a r rank  &fs]/pl  , 


where 


fxl  A <;mallost  integer  m,  such  that  m ^ x 


Ilcnarks 


ii 

[t 

f 


5 


} 

i; 


1.  If  p = 1 (1.0.,  v4ien  S has  scalar  entries) 

then 

a = rank  6[S]  . 

? . \Qien  S Is  a Toeplltz  natrlx,  then  &(S]  =0 
and  therefore  a = 0.  If  S Is  an  arbitrary 
natrlx,  5CS]  nay  be  full  rank  and  then  Ct  = 
N.  Therefore,  a Is  bounded  by  0 < a < N, 
and  the  actual  value  of  a Is  an  index  of 
"distance  from  stationarity"  of  the  matrix  S. 

Kxar.plos 


Lot  T be  a full  Toeplltz  matrix  and  Lj_(Uj^) 
bo  lower  (upper)  Toeplltz  matrices. 

(1)  S = LU  , a < 1 . 

It  can  bo  shevn  (check  for  a 3x3  case)  that, 

JCUJ]  = lu' 

where  1 Is  the  first  column  of  L and  u’  Is 
the  first  row  of  U.  Therefore,  rank  J[LU]  Is 
at  most  p. 

(11)  S = T . a < 2 . 

For  Toeplltz  m.atrlces  6[Tl  = 0,  and  In  J[T] 
only  the  first  row  and  column  may  be  nonzero. 

(Ill)  S = U . L , a < 3 . 

(Iv)  S = T^Tg  , a < 4 ,. 

k 

(v)  s = X)  L ' a < k . 

1=1 

since  the  matrix  &CS]  has  a rank  < a-p  (by 
the  definition  of  a),  wo  can  always  factor  it 
(nor.unicucly)  as 

SCSI  = Pj  E 

with  a sienaturc  matrix  Z and  Pj , P2  being 
NpXOp  (block)  r.atriccs. 

If  S is  a s>Timetric  matrix,  so  is  6[S]. 

In  this  case,  a (nonunique)  svinmctric  decomposi- 
tion can  always  be  found  of  the  form 

6CS]  = P E P’  (1.0. , Pj  = 

where  Z is  again  the  signature  matrix. 

The  question  of  how  to  find  this  decomposi- 
tion and  further  details  about  the  displacement 
rank  a can  be  found  in  [10,11].  For  our  present 
discussion  it  suffices  to  know  that  given  an  ar- 
bltr.ary  covariance  matrix  R of  an  observed  pro- 
cess >■(•),  *c  can  associate  with  it  a number  a 
•and  nr.  Np  X fTi  matrix  P such  tiiat  6[f!)  = PZP'. 
An  Iriporlas.t  example  tliat  further  i 1 Uistra tos  the 
r.ennlng  of  Cl  and  F’  will  he  presented  In  .Sec.  5. 


Wo  can  now  proceed  in  solving  a general  esti- 
mation problem  by  stating  our  assumptions  on  the 

crosscovarlanco  R , 
xy 

B[R  ] = P E P'  (8) 

xy  xy 

where  P Is  a Nn  XOp  matrix, 
xy 

Note  that  if  the  processes  x(«)  and  y(*) 
are  jointly  stationary,  then 

BCR  ] = 0 , BCR)  = 0 , P„„  = 0 . P = 0 . 

xy  *y 

The  motivation  behind  this  assumption  is  that 
in  many  problems  the  signal  x(')  and  the  obser- 
vations y(')  are  connected  by  a linear  relation 
of  the  foi-m 

y^  = Hx^  ''i  " ’'l  '^l 

where  v(*)  is  white  noise  with  unit  intensity. 


'uncorrelated  with  x(’). 

, In  this  case, 

’^i.j  = ’ 

' * ^.j  = *="i^j  ^ ' • 

= HEx^yj  t I . + I • B^_^ 

(10) 

B A j 

[0  1 ><  J 

i 

11 

[l  1 = J 

and  U^y  will  indeed  satisfy  tlie  nssumptlon  (8). 
As  a matter  of  fact, 

R = dlag  (H)  R^  + I 

xy 


and  by  operating  with  &[•]  on  both  sides  we  get 
P = diag  (H)  P . 


4.  The  I/Cvinson-Type  Algorltlin  for 
the  Joint  (x,y)  Process 


Using  the  assumptions  stated  in  the  previous 
section  on  the  covariance  information,  wc  can  now 
give  a set  of  recursions  for  computing 


,m+l  r ,m  T ,,-l.,m’ 
h = 0 ,h  + E M n 

xy  L xyj  m m ' 


.m+1 


m+1 


x“^c 

ra  m 


1 xy  -1 

1 S3  r r 
xy  1,0  0,0 


A°  = I 


(lib) 


B = Cl,0 01 

(11c) 


N , =-•  N’  - c ;.f^c  ' 

m+1  T.\  m n m 


M = M - C’N  C 
in+1  r;  m m v:. 


% = ^0,0 


r 0"1 

0,0  ' 

0 -r 


(12b) 


n 


I 


» i 


I 


1 


Dofluo  a r:  ((t+Dj)  ruul  m = wp  (the  m l:'.o  of  !('"). 
The  illmenslons  of  l>xyi  ^iii  ■ ' Vi  > 'miu 

nxiit,  (m+]))Xp,  (ni  >p)X>’7,  pxp,  clxu,  ii  XU,  i>nd 
p Xu,  rosiK'c'tlvoly.  The  fpiaii 1 1 1 los  E,„,  Cm  linvo 
to  bo  cor.jmtcU  nt  each  step  by 


E = fo  h'"  1 

« L xyj 


■ ["  “J 


*0,0 

0 

*m,0 

'‘0,0 

0 ” 

I 

'‘m+l  ,0 

P E 

• • • 1 

r , 1 

m +1 ,mj 

„,i,'Tn  (ij) 

|_n+l,0  ni  J 


whoro  P|„,  are  t!io  block  row  of  P, 

respectively,  and 


p"  ^ 


By  counting  the  nur.bcr  of  operations  required 
at  the  step  of  the  recursion,  «o  get  (assuming 
p « ni  and  Ignoring  terms  accordingly)  (2n  +3p)IJQ 
liiultipllcatlons . ITndlng  will  tliercforo  I’c- 

qulre  ~(1 .5 +n/p)m“a  nultiplicat  ions  . llie  proof 
of  the  recursions  above  Is  somowl\at  longtliy  and 
Is  given  in  Appendix  A.  It  Is  however, 

important  to  note  that  the  auxiliary  quantities 
A'^,  ft"’  also  obey  the  following  equations. 


.TH 

R A 


(15) 


ni 


0 

pm-1  . 


0 


(16) 


The  first  equation  implies  that 
a"  = [-h^I]’ 

so  that  A™  is  Just  the  optimal  predictor  defined 
in  the  previous  section.  Note  also  that  in  the 
stationary  case 


■n  m 
R U 


in  which  case  II'"  is  the  si.ioolhliig  flltei'Tor  es- 
timating >'0  given  {y^,  1 < 1 •'  n),  oi-  ill  this 
ease  also  Iho  .so-cullo<l  "backward  predictor"  [10^ 
20]. 


Recursive  solutions  of  this  typo  were  devel- 
oped for  the  slationary  ease  by  Lovliisou  [il]  for 
computing  h(‘,")  and  by  Wiggins  and  Robinson 
[15]  for  comijutlng  Indeed,  when  v,e 

take  U = 0,  P =;  0 tiio  equations  (11b),  (lie), 
(12)  reduce  to  the  levinson  algorltlini  and  (11), 
(12)  can  bo  shown  to  bo  equivalent  to  the  equa- 
tions of  iiTggins  and  Robinson.  Thus,  the  station- 
ary ease  is  nicely  imbedded  in  our  framework. 


5.  State  Space  Structure  and 
Cl>andvasokhar-T>'pe  Kquations 

TIio  results  described  so  far  are  quite  gen- 
ernl  and  <Io  not  require  a state  space  structure. 

Wo  shall  now  sliow  how  by  imposing  more  structure 
on  the  covariance  matrices,  the  aiandrasoldiar-typc 
equations  can  bo  derived  from  the  Levinson-typo 
equation  presented  in  the  previous  section. 


Lot  y(")  and  x(*)  be  the  output  and  the 
state  vectors  of  a llneai*  system  driven  by  white 
noise , 1 .e . , 


1+1 


FiXi  + 


(17) 


yi  = H^x^  + 


Eu  u * =5  Q 5 . 

^ i J ' 

Eu  V = Ex  U * =0 
^ 1 j o"j 


Ev  V * =s  1 5 

^°i,j 


In  this  case, 


*^i+l,j  ‘i*i,j  • 


uhere  the  last  term  1«s  zero  for  1 > i.  Thoroforo, 

1 > j . (IS) 

Also,  as  already  noted  earlier, 

r,  , = JI  . + 16.  , . (19) 

i|j  ii,J  ^»J 


In  the  following  discussion,  we  shall  therefore 
assume  that  R,  R^y  obey  assumptions  (IS),  (19) 
which  ai-o  somewhat  weaker  than  the  state-space  as- 
sumption above. 

r.ic  optimal  filter  hjjy(-,0  was  showu  to 
obey  equation  (1).  Tiicrcforc, 


i,s 


1 


xy 

t+1  ,s 


t-1 


r . 
i,s 


0 


Subtractliii;  those  last  equations  and  uslni;  (18), 

VO  got 

t-1  [-  , 

^ Ji  Ct+1,1)  - h (t,i)  +h  (t+l,t)  r 

1^0  L xy  J 1,3  xy  t,s 


t-1  r 

i [" 

i=0  1 


(t+1,1)  - h (t 


"i.s 


xy  • xy  ■ j l,s 

.xy 


nnd  by  comparison  with  (1) 


= (f^  - I - h (t+l,t)  If  ) h (t,l) 

\ t xy  t / xy  ’ 

or 

h (t+l,i)  = (f.  -h  (t+l,t)  h (t,l)  . (20) 

xy  \ t xy  t/  xy  ' 


nnd  Irora  the  Levinson  recursions  (11a),  (lie)  for 

h nnd  U wo  have 
xy 

h (t+l,t)  = 1-  (t,t-l)  + E 

xy  ’ xy  ’ t t \ t/ 

t -1  t t 

= the  last  block  row  of  B , 

so  that 

h^^.(t.i,t)  = ^/t.t-1)  . (23) 

This  recursion  can  be  rewritten  In  another  form 
that  is  easier  to  compare  with  the  usual  Chandra- 
sekhar equations, 

h^y(t+i,t)  = h^y(t.t-i) 

(24) 

The  necessary  algebra  to  derive  this  from  (23)  is 
given  in  Appendix  B. 

In  Appendix  C we  shall  also  prove  that 

*^t  = (^  " ”)  ^-1  • 


Using  the  estimator  equation 


wo  got  a recursive  formula  for  the  estimate, 
t 

x(t+l(t)  = y h (t+1,1)  y 
1=0  ^ 

. = (^t  - \-y<^tl,t)  H^)  Y h^y(t,i)  y^ 

+ h^y(t+l,t) 

= F^x(tlt-l)  +h  (t+l,t) 
t ' xy 

• (y^  -h^y(t+l,t)  H^x(tlt-l))  , (21) 

the  usual  Kalnan-f liter  equation  for  the  state  es- 
timates. Note  that  only  hxy(t+l,t)  is  required 
under  the  state-space  structure  assu.r.ptions ! 

In  the  following,  we  shall  assume  that  F and 
H are  constant.  By  considering  the  defining 
equations  (1),  (2)  of  h(*,‘)  and  h^y(-,'),  it 
is  easy  to  see  that  assumption  (19)  leads  to 

h(t,s)  = IDi  (t,s)  . 

Therefore,  from  the  definitions  (13),  (14)  of 
and  C’lT  it  follows  that 


C 

m 


UK 


n 


Fln.ally,  noto  thnt  (12),  (24),  (25)  pi-ovldc  a com- 
ploto  set  of  recursions  for  computing  h.{y (t+1  ,t) , 
which  are  of  the  CliandrasoUhar-typo . Indeed,  a 
comparison  of  our  rc.sults  to  those  presented  in 
[IGI  shows  that  equations  (24),  (25),  (12a),  (12b) 
arc  precisely  equations  (16),  (13'),  (15),  (14)  of 
[16]  if  the  following  chango  of  notation  is  made. 

h (t+l,t)  <->  K (t)  (or  h (t+l,t)  N <->  K 1 , 
xy  ' g I xy  ’ ' t t/  ’ 

N.  <->  R^(t)  , M,.  <->  -if  , <->  y^  . 

V V X t t 

We  see  therefore  that  the  Qiandrase’idiar  equa- 
tions arc  naturally  induced  by  the  Levinson-typo 
recursions  when  the  coi-ariancc  matrices  have  a 
special  structure.  Furthermore,  the  parameter  C<o 
which  appears  in  the  Cliandraseldiar  algorithm  car. 
now  be  shown  to  have  a meaning  in  terras  of  Q (or 
S) , the  displacement  rank  of  R (sec  also  fS,9]). 

To  see  this,  let  us  recall  that  for  a tlmc- 


Invariant 

state-space  model 

J 

i > 

3 

*-i.j  “ W’ 

i = 

3 

(26) 

i < 

3 

where 

and 

''i.i  = 

+ Q . 

(27) 

so  that , 

^i+1 

- q„,  ■ 

,1-1  - ’^.1' 

- n„)(rS 

'll’  . 

(22) 


WrltliiH;  tills  In  matrix  fomi  i;lvo3 
= 0 (itj  - i(p)  O' 


II 

OA 


nn  oxtonclocl 
, obsorvnblllty 
matrix . 


let  us  take  for  simplicity  the  scalar  case  (l.o., 

P = li  r^  j scalar).  Assunlnc  that  the  system  is 
obsei-a’ablo , O will  be  a full  rank  matrix,  and 
therefore 


a = rank  6H  = rank  (rt^  - it^) 


Now  the  Chandrasekliar  equations  involve  a parame- 
ter Oq  defined  as 

Oq  = rank  - it^)  - FUqII’ + D'SwqF'J 

and  since  the  second  term  In  the  brackets  is  of 
rank  1, 

n-i<ao<a  + l . 


Finally,  wo  should  note  tJiat  Ihoso  d Isci-oto- 
tlmo  results  have  close  continuous -time  atmlo;;s 
presented  In  [13]  for  the  i;enoral  preblco  of  solv- 
liip  some  Intci-jral  equations,  and  In  [12]  deallni;  ' 
specifically  with  Uio  estimation  probleii.  In  fact, 
it  was  these  results  that  provided  the  Imiiodlate 
motivation  for  the  dlscroto-timo  analysis  pre- 
sented hoi’o , 


Appendix  A 

Proof  of  the  Lovlnson-Typo  Algorithm 
for  the  Joint  (x,y)  Process 

The  defining  cciuation  for  h"^  was  given  ns 

xy 


h'"  r"-"  = ...  1 

xy  L '"i®  ra,m-lj 


Therefore , 
. m 1 

° ’■.yj  " 


0 ...  0 ...  0 I r r 


Fo  h"  1)  : 

L xyj<;  . m-1 


0,0 *‘0,m 


pm-lj.pm-1 


6.  Conclusions 

Vo  have  shown  how  recursive  solutions  can  be 
obt-ained  for  the  optimal  predictor  witlc  covariance 
(or  input-output)  data,  whether  or  not  statc-spaco 
models  are  available.  Tlic  complexity  of  these  al- 
gorithms depends  on  a measure  Q.  of  the  "distanco" 
fror.  statlonnrity  of  the  signal  and  observed  pro- 
cesses . 

A similar  ajiproach  makes  It  possible  to  dci’ivo 
a recursive  solution  for  the  optimal  (filter) 
smoother  H.vy  The  details  will  not  be  presented 
hero  (see  [10]  for  a partial  treatment),  but  It  is 
important  to  note  that  Hjjy  can  bo  computed  in 
0{CC:~)  operations,  instead  of  O(N^)  required  for 
a direct  solution  of  equation  (4). 


... 

J + \l 

V ,0  ...  o| 

L ni,0 

n,m-lj  L 

f J 

/ m m-l„> 

1 r m-ll 

+ (h  P L, 
\ xy  > 

|[0.P  J 

F^m+1,0  '■  n+l,isj  L m+l’**j!i  J 


Fa,  ,0  . 

/ m m-l„\r 

m-ll 

. . 0 + ( h p ) 

O.P 

L 1 

J V xy  /L 

J 

r 1 ° 

^ = [°‘’xy]  : 


In  the  special  case  where  the  processes  are 
known  to  come  from  a constant-pai'amctcr  statc- 
spacc  r.ioclcl  , the  distance  from  stationarity  O, 
coincides  with  a parameter  describing  the  computa- 
tional reductions  obtainable  by  using  tlic  previ- 
ously known  [IG]  Oiandrasoidiar  equations.  '.lore- 
over,  those  general  recursions  reduce  naturally  to 
the  Cnandrascldiar  equations  in  this  and  actually 
also  in  some  more  general  cases.  Note,  for  exam- 
ple, that  wo  made  no  assumption  on  C)  and  our 
derivation  holds  when  it  is  time  varying.  Note 
ftl '.o  that  our  approach  leads  to  a derivation  of 
tlio  Chandrasekhar  equations  that  does  not  mention 
the  l;icc.ati  ctiuatlons,  v.iiich  was  at  the  heart  of 
the  original  derivation  [10 1 . Aelually,  the  gen- 
eral tlnc-vnrying  Rlccntl  equation  can  also  he  Im- 
bediled  In  tlic  framework  presented  hero  (see  [17]), 


This  can  be  rewritten  as 


r a , fl  0 . . . O" 

Fo  h"  1 r"  = h^'V  - F ' m-1'  (A2) 

L xyj  xy  P 


r 0 

0,0 


m = L xyj  . L “+1-  ■'  J 


Define  auxiliary  quantities  A"’ , D"  as 


(a:i) 


•vliovo  Iht.  liilit  bloc!;  row  tif  A'*'  Is  tho  Idoiitlty 


I 0 


, n n -1  0 , 

U I!  M = . ,n-l 

m . H 


Wito  t!:at  froT  (A2)  nnd  (A4)  It  follows  that 

= [o  h'"  1 + E J,rV’  (A5) 

.1  xy  -1 

and  tro\r\  the  deflnln:'  equation,  = ^1, 0^0,0* 

t'slni;  a sinllar  ;.pproacii,  tho  rccuri>lons  for 
are  derived. 


-■'"j  I-:., 


0,n+l 


M : p"j:i 


“ 0 * p"  ’ 


M = M - F’N  C , 
m+1  m m m m 0 


r 0“> 

0.0 


0 £ 


. (All) 


To  vci-lfy  (A8),  (AO),  wo  preniultlply  tlic  equations 
by  and  clioeU  Uiat  tho  rli;ht-hand  side  of 

those  equations  satisfy  tho  definin';  equations 

(a;j),  (ai). 

Finally,  It  remains  to  show  that  = F^, . 

Tlie  proof  is  lengthy  and  shall  be  omitted.  It  can 
be  found  in  [10,11). 


Appendix  B 

An  Alternative  form  of  the  Kecurslon  for  h (t+1  t) 

xy 

i.,y(t>i.t)  = (\/t,t-i)  - v't''^;-i"’-"t-i) 

• (\-i  - '"^t-A-i";-i") 


'"o.o  ” 


P.-^^  M-"=  : 

0 " °-- 


k 

m m 


- h <t,t-i)  HE.  e;  H’ 

xy  t-1  t-1  t-1 


^ vC(^t-i»'Vi”^-i " -^-i) 


k ^ fr  r _ 1 b"' 

m ~ m+1,0  m+l,mj 

Combinlnp:  the  defining  equations  of  A™,  B*^  and 
(A6),  (A7)  gives 

1 r°  1 

= - m"  F'  , A-  = I . (A8) 

.n  . n m 0 

A 0 


N , = N - (k  - [0,p  1 M ) M"^f 
m+1  r.i  V n ’ m m ) n 


N ' 0 M F*  , = r„  „ 

n in  n r.  ’ 0 0,0 


= h^y(t,t-l) 


In  Appendix  C,  wo  shall  show  that 


E^  = (f  - h (t,t-l)  ll)  E,  , 
t ( xy  /t-1 


there fore 


h (t+l,t)  N.  =h  (t,t-l)  , - FE,  ,E’  II'. 

xy  t xy  ’ t-1  t-1  t-1  v-1 


C k - ro,p  1 M 
n ~ m n n 


1 -1 
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ApiK*ndix  C 


The  Kecur.slon  of  K 


An  a firj.t  stoii,  v.i  ..Ijall  i);o.v 


I. 


h (t,s)  = (?  ~ h ll)  h 

xy  ’ \ y-y  / xy 


0 < s < t-2 


(Cl) 


KroT  tho  dofinlnjc  eciuatlou  (1)  for  h , wo  know 


that 


xy 


("Jv  - 

■fe 


where 


Uslnir 


gives 


O.t-1 


L t-i.t-i. 


= Fi-’'y 

t,s  t-l,s 


■ 

[“ oaJ 

where  Is  whatever  is  necessary  to  satisfy  tho 

last  equation.  Hence, 

-xy  - - ‘>[-xy’“] 

+ \.y(t  pt-l)  |^-h*"^lj  (C2) 

Usir.3  h(t,s)  = }fh  Ct,s),  v.e  can  rcwTlto  (C2)  as 
xy 

= (f  - h (t,t-l)  ol 

xy  \ xy  ' / L J 

+ [o  ...  0,h^y(t,t-l)] 

which  proves  (Cl) . 

Let  us  rcwTltc  the  definition  (13)  of  as 


0,0 


t-2 

P L 


t-1,0 
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(C3) 


It  is  also  true  that 


Using  (Cl),  (Ci),  and  (C5) , wo  can  rov/rlto  (C3)  as 
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>1  = [f  - h (t,t-l)  nl 
t L xy  ’ J xy  I 
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- ('  - -xy<*’—>  ”)  ^-1  • 
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